<html>

<head>
<meta http-equiv="Content-Language" content="en-us">
<meta name="GENERATOR" content="Microsoft FrontPage 5.0">
<meta name="ProgId" content="FrontPage.Editor.Document">
<meta http-equiv="Content-Type" content="text/html; charset=windows-1252">
<title>Twiss Loop</title>
</head>

<body>

<font SIZE="2">
<address>
</address>
<p align="center"><font SIZE="2" COLOR="#0000ff"><b>SUBROUTINE</b></font> 
compute_twiss(r,my_state,filename,pos,del,no,thinlens,name) of PTC</p>
<p align="center">in Sra_fitting.f90</p>
<address align="center">
&nbsp;</address>
<address>
  STATE=my_state+delta0
</address>
<address>
  closed=zero</address>
<address>
  closed(5)=del</address>
<hr>
<address>
  </font><font SIZE="2" COLOR="#0000ff"><b>CALL</b></font><font SIZE="2"> 
  FIND_ORBIT(R,CLOSED,pos,STATE,c_1d_5)&nbsp;&nbsp;
</address>
<hr>
<address>
  </font><font SIZE="2" COLOR="#0000ff"><b>write</b></font><font SIZE="2">(6,*)
  </font><font SIZE="2" COLOR="#800000">&quot;closed orbit &quot;</address>
</font><font SIZE="2">
<address>
  </font><font SIZE="2" COLOR="#0000ff"><b>write</b></font><font SIZE="2">(6,*) 
  CLOSED</address>
<address>
  </font><font SIZE="2" COLOR="#0000ff"><b>CALL</b></font><font SIZE="2"> 
  INIT(STATE,no,0)</address>
<address>
  </font><font SIZE="2" COLOR="#0000ff"><b>CALL</b></font><font SIZE="2"> 
  ALLOC(Y); </font><font SIZE="2" COLOR="#0000ff"><b>CALL</b></font><font SIZE="2"> 
  ALLOC(NORM)</address>
<address>
  </font><font SIZE="2" COLOR="#0000ff"><b>call</b></font><font SIZE="2"> 
  alloc(id,a_f,a_l,a_nl,DR,R_TE,cs_te)</font></address>
<hr><font SIZE="2">
<address>
  id=1</address>
<address>
  Y=CLOSED+id</address>
<address>
</address>
<address>
  </font><font SIZE="2" COLOR="#0000ff"><b>CALL</b></font><font SIZE="2"> 
  TRACK(R,Y,pos,STATE)</address>
<address>
  id=y</address>
<address>
  NORM=ID&nbsp;&nbsp;&nbsp; <span style="font-style: normal"><b>
  <font color="#FF00FF">! Normalize the one-turn map from position pos back to 
  pos</font></b></span></address>
<hr></font>
<address>
  <font SIZE="2" COLOR="#0000ff"><b>write</b></font><font SIZE="2">(mf,*) </font>
  <font SIZE="2" COLOR="#800000">&quot; fractional tunes &quot;</font><font SIZE="2">,norm%tune(1:2)</address>
<address>
  p=&gt;r%start</address>
<address>
  s=zero</address>
<address>
  </font><font SIZE="2" COLOR="#0000ff"><b>do</b></font><font SIZE="2"> i=1,pos-1&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;
  <span style="font-style: normal"><b><font color="#FF00FF">! </font></b></span>
  <font color="#FF00FF"><b><span style="font-style: normal">Adds up the lengths 
  from s=0 to s(pos)</span></b></font></address>
<address>
  s=s+p%mag%p%ld</address>
<address>
  p=&gt;p%next</address>
<address>
  </font><font SIZE="2" COLOR="#0000ff"><b>enddo&nbsp; </b></font>
  <font SIZE="2"><span style="font-style: normal"><b><font color="#FF00FF">!
  </font></b></span><font color="#FF00FF"><b><span style="font-style: normal">
  moves fibre pointer</span></b></font></font></address>
<hr><font SIZE="2">
<address>
  px=zero&nbsp; <span style="font-style: normal"><b><font color="#FF00FF">!&nbsp;
  </font></b></span><font color="#FF00FF"><b><span style="font-style: normal">
  Phase advance is set to zero</span></b></font></address>
<address>
  py=zero</address>
<address>
</address>
<address>
  </font><font SIZE="4" COLOR="#FF0000"><b>call 
  factor(norm%a_t,a_f,a_l,a_nl,DR=dr,r_te=r_te,cs_te=cs_te,te_done=te_done)&nbsp;</b> </font><font SIZE="2">
  <span style="font-style: normal"><b><font color="#FF00FF">!&nbsp; Calls the 
  routine to factor the map and put it in a canonical form</font></b></span></address>
<address>
  Y=CLOSED+norm%a_t</address>
<hr></font><font SIZE="2" COLOR="#008000">
<address>
  !&nbsp;
</address>
</font><font SIZE="2">
<address>
  bx=(norm%a_t%v(1).sub.'1')**2+(norm%a_t%v(1).sub.'01')**2
  <span style="font-style: normal"><b><font color="#FF00FF">!&nbsp; </font></b>
  </span><font color="#FF00FF"><b><span style="font-style: normal">beta_x</span></b></font></address>
<address>
  by=(norm%a_t%v(3).sub.'001')**2+(norm%a_t%v(3).sub.'0001')**2
  <span style="font-style: normal"><b><font color="#FF00FF">!&nbsp; </font></b>
  </span><font color="#FF00FF"><b><span style="font-style: normal">beta_y</span></b></font></address>
<address>
  ax=-(norm%a_t%v(1).sub.'1')*(norm%a_t%v(2).sub.'1')-(norm%a_t%v(1).sub.'01')*(norm%a_t%v(2).sub.'01')&nbsp;&nbsp;
  <span style="font-style: normal"><b><font color="#FF00FF">! </font></b></span>
  <font color="#FF00FF"><b><span style="font-style: normal">alpha_x</span></b></font></address>
<address>
  ay=-(norm%a_t%v(3).sub.'001')*(norm%a_t%v(4).sub.'001')-(norm%a_t%v(3).sub.'0001')*(norm%a_t%v(4).sub.'0001')&nbsp; 
  !<span style="font-style: normal"><b><font color="#FF00FF"> </font></b></span>
  <font color="#FF00FF"><b><span style="font-style: normal">alpha_y</span></b></font></address>
<address>
  ex=a_f%v(1).sub.'00001'&nbsp; <span style="font-style: normal"><b>
  <font color="#FF00FF">&nbsp;! </font></b></span><font color="#FF00FF"><b>
  <span style="font-style: normal">eta_x</span></b></font></address>
<address>
  exp=a_f%v(2).sub.'00001' <span style="font-style: normal"><b>
  <font color="#FF00FF">&nbsp;&nbsp; ! </font></b></span><font color="#FF00FF">
  <b><span style="font-style: normal">eta_y</span></b></font></address>
<address>
  phi=(r_te%v(1).sub.'1')-one <span style="font-style: normal"><b>
  <font color="#FF00FF">! </font></b></span><font size="2" color="#FF00FF"><b><span style="font-style: normal">
  Teng-Edwards cosine or hyperbolic cosine</span></b></font></address>
<address>
  <font SIZE="2" COLOR="#0000ff"><b>if</b></font>(<font SIZE="2" COLOR="#0000ff">abs</font>(phi)&lt;epsflo) 
  phi=zero</address>
<address>
  <font SIZE="2" COLOR="#0000ff"><b>if</b></font>(COSLIKE)
  <font SIZE="2" COLOR="#0000ff"><b>then</address>
</b></font>
<address>
  typec=<font SIZE="2" COLOR="#800000">' = COS-1 '</address>
</font>
<address>
  <font SIZE="2" COLOR="#0000ff"><b>else</address>
</b></font>
<address>
  typec=<font SIZE="2" COLOR="#800000">' = COSH-1'</address>
</font>
<address>
  <font SIZE="2" COLOR="#0000ff"><b>endif</address>
</b></font>
<address>
  <font SIZE="2" COLOR="#0000ff"><b>write</b></font>(mf,<font SIZE="2" COLOR="#800000">'(6x,a4,6x,(6x,a1,6x),2(4x,a5,4x),2(4x,a6,3x),2(4x,a6,3x),2(3x,a7,3x),a25)'</font>) 
  &amp;</address>
<address>
  <font SIZE="2" COLOR="#800000">&quot;name&quot;</font>, <font SIZE="2" COLOR="#800000">
  &quot;s&quot;</font>, <font SIZE="2" COLOR="#800000">&quot;betax&quot;</font>,<font SIZE="2" COLOR="#800000">&quot;betay&quot;</font>,<font SIZE="2" COLOR="#800000">&quot;alphax&quot;</font>,<font SIZE="2" COLOR="#800000">&quot;alphay&quot;</font>,<font SIZE="2" COLOR="#800000">&quot;eta_x 
  &quot;</font>,<font SIZE="2" COLOR="#800000">&quot;etap_x&quot;</font>,<font SIZE="2" COLOR="#800000">&quot;Phase 
  x&quot;</font>,<font SIZE="2" COLOR="#800000">&quot;Phase y&quot;</font>,<font SIZE="2" COLOR="#800000">&quot;Teng-Edwards 
  Cos or Cosh &quot;</font></address>
<address>
  <font SIZE="2" COLOR="#0000ff"><b>write</b></font>(mf,<font SIZE="2" COLOR="#800000">'(a16,10(1x,g12.5),1x,a9)'</font>) 
  p%mag%name(1:16),s,bx,by,ax,ay,ex,exp,px,py,phi,typec</address>
</font>
<hr>
<address>
  <font SIZE="2">
</address>
<address>
</address>
<address>
  </font><font SIZE="2" COLOR="#0000ff"><b>do</b></font><font SIZE="2"> i=pos,pos+r%n-1&nbsp;&nbsp;
  <span style="font-style: normal"><b><font color="#FF00FF">! Do looping over 
  the fibre</font></b></span></address>
<hr>
<address>
  &nbsp;&nbsp;&nbsp;
  t=&gt;p%t1</address>
<address>
  </font><font SIZE="2" COLOR="#0000ff"><b>do</b></font><font SIZE="2"> </font>
  <font SIZE="2" COLOR="#0000ff"><b>while</b></font><font SIZE="2">(.not.</font><font SIZE="2" COLOR="#0000ff">associated</font><font SIZE="2">(t,p%t2%next))&nbsp;
  <span style="font-style: normal"><b><font color="#FF00FF">! Do looping over 
  integration nodes</font></b></span></font></address>
<address>
  <font SIZE="2" COLOR="#0000ff"><b>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; if</b></font><font SIZE="2">(thinlens)
  </font><font SIZE="2" COLOR="#0000ff"><b>then</address>
</b></font><font SIZE="2">
<address>
  </font><font SIZE="2" COLOR="#0000ff"><b>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; call</b></font><font SIZE="2"> 
  track_probe_x(y,STATE, node1=t,node2=t%next)&nbsp;&nbsp;
  <span style="font-style: normal"><b><font color="#FF00FF">! Do looping over 
  the fibre; sometimes gets the integer part of the tune wrong (too strong 
  magnet)</font></b></span></font></address>
<address>
  <font SIZE="2" COLOR="#0000ff"><b>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; else</address>
</b></font><font SIZE="2">
<address>
  </font><font SIZE="2" COLOR="#0000ff"><b>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; call</b></font><font SIZE="2"> 
  track_probe_x(y,STATE, node1=t,node2=p%t2%next)
  <span style="font-style: normal"><b><font color="#FF00FF">! Do looping over 
  the thin lens</font></b></span></font></address>
<address>
  <font SIZE="2" COLOR="#0000ff"><b>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; endif</address>
<hr></b></font><font SIZE="2">
<address>
  closed=y</address>
<address>
  norm%a_t=y</address>
<address>
  </font><font COLOR="#FF0000"><b>call</b></font><b><font color="#FF0000"> 
  factor(norm%a_t,a_f,a_l,a_nl,DR=dr,r_te=r_te,cs_te=cs_te,te_done=te_done)</font><font SIZE="2">&nbsp;&nbsp;</font></b><font SIZE="2">
  <span style="font-style: normal"><b><font color="#FF00FF">!&nbsp; Calls the 
  routine to factor the map and compute the phase advance</font></b></span></address>
<address>
&nbsp;px=px+</font><font SIZE="2" COLOR="#0000ff">asin</font><font SIZE="2">(dr%v(1).sub.'01')/twopi
  <span style="font-style: normal"><b><font color="#FF00FF">! </font></b></span>
  <font color="#FF00FF"><b><span style="font-style: normal">phase x</span></b></font></address>
<address>
&nbsp;py=py+</font><font SIZE="2" COLOR="#0000ff">asin</font><font SIZE="2">(dr%v(3).sub.'0001')/twopi
  <span style="font-style: normal"><b><font color="#FF00FF">! </font></b></span>
  <font color="#FF00FF"><b><span style="font-style: normal">phasey</span></b></font></address>
<hr>
<address>
  Y=CLOSED+norm%a_t</address>
<address>
  </font><font SIZE="2" COLOR="#0000ff"><b>if</b></font><font SIZE="2">(thinlens)
  </font><font SIZE="2" COLOR="#0000ff"><b>then</address>
</b></font><font SIZE="2">
<address>
  t=&gt;t%next</address>
<address>
  </font><font SIZE="2" COLOR="#0000ff"><b>else</address>
</b></font><font SIZE="2">
<address>
  t=&gt;p%t2%next</address>
<address>
  </font><font SIZE="2" COLOR="#0000ff"><b>endif</address>
</b></font><font SIZE="2">
<address>
  </font><font SIZE="2" COLOR="#0000ff"><b>enddo</b></font><font SIZE="2">
</address>
<address>
  s=s+p%mag%p%ld</address>
<address>
  p=&gt;p%next</address>
<address>
  bx=(norm%a_t%v(1).sub.'1')**2+(norm%a_t%v(1).sub.'01')**2
  <span style="font-style: normal"><b><font color="#FF00FF">!&nbsp; </font></b>
  </span><font color="#FF00FF"><b><span style="font-style: normal">beta_x</span></b></font></address>
<address>
  by=(norm%a_t%v(3).sub.'001')**2+(norm%a_t%v(3).sub.'0001')**2
  <span style="font-style: normal"><b><font color="#FF00FF">!&nbsp; </font></b>
  </span><font color="#FF00FF"><b><span style="font-style: normal">beta_y</span></b></font></address>
<address>
  ax=-(norm%a_t%v(1).sub.'1')*(norm%a_t%v(2).sub.'1')-(norm%a_t%v(1).sub.'01')*(norm%a_t%v(2).sub.'01')&nbsp;&nbsp;
  <span style="font-style: normal"><b><font color="#FF00FF">! </font></b></span>
  <font color="#FF00FF"><b><span style="font-style: normal">alpha_x</span></b></font></address>
<address>
  ay=-(norm%a_t%v(3).sub.'001')*(norm%a_t%v(4).sub.'001')-(norm%a_t%v(3).sub.'0001')*(norm%a_t%v(4).sub.'0001')&nbsp; 
  !<span style="font-style: normal"><b><font color="#FF00FF"> </font></b></span>
  <font color="#FF00FF"><b><span style="font-style: normal">alpha_y</span></b></font></address>
<address>
  ex=a_f%v(1).sub.'00001'&nbsp; <span style="font-style: normal"><b>
  <font color="#FF00FF">&nbsp;! </font></b></span><font color="#FF00FF"><b>
  <span style="font-style: normal">eta_x</span></b></font></address>
<address>
  exp=a_f%v(2).sub.'00001' <span style="font-style: normal"><b>
  <font color="#FF00FF">&nbsp;&nbsp; ! </font></b></span><font color="#FF00FF">
  <b><span style="font-style: normal">eta_y</span></b></font></address>
<address>
  phi=(r_te%v(1).sub.'1')-one <font color="#008000">&nbsp;</font><span style="font-style: normal"><b><font color="#FF00FF">! </font></b></span><font size="2" color="#FF00FF"><b><span style="font-style: normal">
  Teng-Edwards cosine or hyperbolic cosine</span></b></font></address>
<address>
  <font SIZE="2" COLOR="#0000ff"><b>if</b></font>(<font SIZE="2" COLOR="#0000ff">abs</font>(phi)&lt;epsflo) 
  phi=zero</address>
<address>
  <font SIZE="2" COLOR="#0000ff"><b>if</b></font>(COSLIKE)
  <font SIZE="2" COLOR="#0000ff"><b>then</address>
</b></font>
<address>
  typec=<font SIZE="2" COLOR="#800000">' = COS-1 '</address>
</font>
<address>
  <font SIZE="2" COLOR="#0000ff"><b>else</address>
</b></font>
<address>
  typec=<font SIZE="2" COLOR="#800000">' = COSH-1'</address>
</font>
<address>
  <font SIZE="2" COLOR="#0000ff"><b>endif</address>
</b></font>
<address>
</font>
  <font SIZE="2" COLOR="#0000ff"><b>if</b></font><font SIZE="2">(printres)
  </font><font SIZE="2" COLOR="#0000ff"><b>write</b></font><font SIZE="2">(mf,</font><font SIZE="2" COLOR="#800000">'(a16,10(1x,g12.5),1x,a9)'</font><font SIZE="2">) 
  p%mag%name(1:16),s,bx,by,ax,ay,ex,exp,px,py,phi,typec</address>
</font>
<address>
  <font SIZE="2" COLOR="#0000ff"><b>enddo </b></font><font SIZE="2">
  <span style="font-style: normal"><b><font color="#FF00FF">! End of Do looping 
  over the fibre</font></b></span></font></address>
<font SIZE="2"><hr>
<address>
  </font>
</address>

</body>

</html>